Locality properties of the free energy fluxes in gyrokinetic turbulence 
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The nature of the nonlinear interactions in gyrokinetic (GK) turbulence, driven by an ion-temperature gra- 
dient instability, is investigated using numerical simulations of single ion species plasma in three-dimensional 
flux tube geometry. To account for the level of separation existing between scales involved in an energetic inter- 
action, the degree of locality of the free energy scale flux is analyzed employing Kraichnan's infrared (IR) and 
ultraviolet locality functions. Due to the nontrivial dissipative nature of GK turbulence, an asymptotic level for 
the locality exponents, indicative of a universal dynamical regime for GK's, is not recovered and an accentuated 
non-local behavior of the IR interactions is found instead, in spite of the local energy cascade observed. 
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Introduction. — The main characteristic of physical flows is 
given by the existence of couplings between different scales of 
motion (£), described mathematically by nonlinear terms in re- 
spect to the dynamical quantities. As result of these couplings, 
the flow will develop a turbulent state for a sufficiently large 
interval of excited scales. While the range of scales available 
to the flow depends on the boundary conditions and on the 
nature of the flow itself (mathematically described by linear 
terms and externally given parameters), the redistribution of 
information among different scales is only due to the nonlin- 
ear terms. As such, the redistribution mechanism is expected 
to have a universal behavior for intervals of scales for which 
the linear terms are negligible, i.e. the inertial zone. This 
picture stands at the basis of the study of turbulence for hy- 
drodynamical and electrically conductive fluids. In principle, 
gyrokinetic (GK) turbulence makes no exception to this pic- 
ture (H 0], albeit with a series of complications due to the 
nature of the linear terms. 

The gyrokinetic formalism is pertinent for the study of 
multi-species plasmas in the presence of strong magnetic 
guide fields B- By eliminating exactly the gyration phase 
of charged particles around the magnetic field lines |4l|5|], the 
dynamical space can be reduced from six dimensions to five. 
From the start, it is seen that the constraint imposed by the 
magnetic guide field on the charged flow creates an anisotropy 
in the system {£ — {£±, In }). The scaling in the two directions 
are linked as result of causality, a hypothesis known as criti- 
cal balance, recently used in the scaling of plasma turbulence 
10,0]. The current work concentrates on the analysis of the 
perpendicular spatial structures (£±) of the fluctuations. 

To understand the dynamics introduced by the nonlinear 
term, the scale redistribution of free-energy (a GK ideal in- 
variant, i.e. a global quantity that remains constant in time 
in the absence of source and sink effects) is usually investi- 
gated. Different works reported that the exchange of energy 
takes place between closest neighbor dyadic structures [8-10]. 
However, although the energy exchanges are local, the ques- 
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tion regarding the locality of the interactions was never ad- 
dressed, i.e. the fact that the local exchanges of energy might 
be generated by the interaction of highly separated (non-local) 
scales. In this work, we describe a quantitative way of as- 
serting the degree of locality for GK turbulence. The idea of 
locality can be seen as the disparity of scales contributing to 
a nonlinear interaction. For an given energy flux through a 
scale, the degree to which each scale contributes to the men- 
tioned flux represents a useful assertion of locality of interac- 
tions, UHtHl]. For the interaction to be local, the contribution 
of highly separated scales should be small and decrease fast 
with the increase in separation. 

For fluid turbulence, the separation of the forced and dis- 
sipative scale ranges leads to the existence of a natural iner- 
tial range that possesses an unique locality exponent, which 
asymptotes dynamically to a 4/3 value II 1311 . However, in the 
case of GK turbulence and for fusion plasma configurations 
in particular, each instability that generates a time unstable 
mode (energy source) is also accompanied by an ensemble of 
time stable modes (energy sinks). These stable modes, known 
as linear damped eigenmodes, are nonlinearly coupled to the 
unstable modes and are responsible for an additional energy 
dissipation route ifTill . This dissipation mechanism acts at a 
scale comparable to the forcing and as such, does not require 
the existence of a classical nonlinear cascade. Thus, the exis- 
tence of a nonlinear cascade process is made possible only if 
the energy injected by the unstable modes is greater than the 
energy dissipated locally by the linear damped eigenmodes. 
Moreover, since the dissipation tends to permeate strongly the 
forced range and take over from the force as the dominant 
effect for smaller scales, characteristic time wise, the nonlin- 
ear cascade process occurs inside a dissipation range Si. 
This represents a novel and nontrivial complication compared 
to fluid turbulence, which leads to a dampening of contribu- 
tions the to energy flux. As result, a unique locality exponent 
can not be found, in spite of the local energy exchanges. 

Gyrokinetic simulations. — In the present work, numerical 
solutions of the nonlinear gyrokinetic equations in flux-tube 
(s — a) geometry IU7I1 are analyzed using a field-aligned coor- 
dinate system (x, y, z, v\\ , fi) with (128, 64, 16, 32, 8) points in 
each direction, respectively. The solutions are obtained by the 
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use of the GENE code, [18] for ion-temperature gradient (ITG) 
driven GK turbulence with physical parameters correspond- 
ing to the Cyclone Base Case Ill9ll . For a better understanding 
of the nonlinear dynamics involved, the analysis is limited to 
the simple scenario of electrostatic fluctuations generated by 
a single ion species (the species index will be omitted) and 
adiabatic electrons, in the context of a large aspect-ratio cir- 
cular cross-section equilibrium model, for which the equilib- 
rium magnetic field is Bq (in units of the magnetic field value 
on the magnetic axis). For details see Ref. [151. 

Considering that the total ion distribution function F is 
split into an appropriately normalized Maxwellian part Fq = 
7r -3/2 e -(«n+MS ) an( j a p ermr bed part /, the non-adiabatic 
contribution of the ion distribution function is given as h = 
f + (Zcf>i/To)Fo, where cf>i is the gyro-averaged electrostatic 
potential, Tq is the ion temperature (normalized to the electron 
temperature) and Z is the electric charge. The time (t) evolu- 
tion equation for the perturbed distribution function reads, 

|[ = G[f] + L c [f] + L„ [/] + N[f, /] + D[f] , (1) 



In the field aligned coordinates, the x label refers to the flux 
surface, while y identifies different field lines laying on the 
same flux surface. Each scales of length (£±) can now be 
easily identified by the norm (k ~ £~7) of the wave-vector 
(k = )s.±) based in the k x , k y space. 

As for any quadratic quantity, the free energy spectral den- 
sity can be considered {£ — ^ k £ k ), for which the balance 
equation reads, 
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where in the rhs of Eq. the terms A k = {^ k , £ k , 2? k } are 
found as, A k = j d8 ■p I /i k yl~ k using the spectral form of 
the evolution equation for the perturbed distribution function 
Eq. CQ, with A k = {G k ,L k = i k + L^,D k }. While in 
Eq. © the linear quantities £, T> and Q are defined involving 
only k local modes and their complex conjugates, for the term 
generated by the nonlinear product T k different modes enter 
in the definition, 



rk = EE rk ' p ' q5k +p+ q 
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Fo^ is the 



where G[f] = - uj n + (v^ + fiB Q - w r 
contribution of the normalized background density (pj n — 
—Rd\ogno/dx) and temperature (ujt = —Rd\ogTo/dx) 
gradients, with R being the the major radius. It represents 
the driving mechanism for GK turbulence and it is respon- 
sible for the injection of free energy into the system. The 
second term appears due to magnetic curvature, Lc[f] — 
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the s — a geometry (a = 0) employed here. The third 
terms contains the parallel dynamics involving magnetic trap- 
ping and linear Landau damping/pumping effects, L\\[f\ = 
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nE>o , h , where vt = y/2To/m is the ion 

thermal velocity, m is the ion mass and the Poisson bracket 
structure is defined as: [f,g] a b = 9 a fdbg — dbfd a g. The 
second to last term represents the nonlinear term, N[f, f] = 
— [4>i, h\ xy , while the last term in Eq. ((TJ contains the dissipa- 
tive effects. To avoid the computational cost of collisional op- 
erators, the dissipation terms have a simple hyper-diffusivity 
form D[f] = —(a z d™ + a V]] 9" )/, where n = 4 and the a's 
parameters are adapted to the problem at hand 1191, 1 14, 1 1 511 . A 
2/3 dealiasing method is used for the nonlinear term. 

To close the equation system, the self-consistent electro- 
static potential is obtained by solving the gyrokinetic Pois- 
son equation, ^^[1 — r o (6)]0i = ZnB no j Jo(A)/du||d/x, 
where the Bessel function J , r (6) = e b I (b) and the mod- 
ified Bessel function Io have the arguments defined as, A = 
^/JlB^k±v T /n and b = k\v^/(2n 2 ), with fi = ZB„/(mc) 
and k± being the perpendicular wavenumber. 

Free energy balance. — In this formulation, the global 
free energy contained in the system is defined as, 8 = 
\ J da;dyd0 ^ hf , where d0 = (TrBano)dzdv\\d(j,. To an- 
alyzed the excitation degree of perpendicular turbulent scales, 
an integral over the d0 infinitesimal element and a Fourier 
decomposition of the remaining (x, y) space are performed. 



The discreet Delta Dirac <5k+ P + q selects only interactions that 
occur between a triad of modes which obey the resonance con- 
dition, k+p + q = 0. The transfer that takes place for a single 
triad, known as the triad transfer, is defined as, 
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where the symmetry in modes q and p is written explicitly 
ll20ll . At the triad level, the free-energy conservation by the 
nonlinear interaction can be written as, T k,p ' q + T p q,k + 

-y-q,k,p _ q 

Although the triad transfers contain the complete physi- 
cal information related to the energetic coupling of scales, 
the sheer number of transfers involved makes them unus- 
able in any direct manner. To ease our work, we decompose 
the spectral space into a series of structures and analyze the 
transfers that occur among them. The structures boundaries 
sk = (kx-i, fcif ] are typically given as a power law in terms 
of the wavenumber k, k^ = k x 2 I - K ~ 1 ^ A . The filtered ion 
distribution function h K and filtered electric potential <j) K are 
found to be: 



{/ l ,<^ 1 } K (k) = 



{h, <M(k), k e s K 
0, k 4 s K 



(5) 



Depending on the selection of ko and A, we obtain a N num- 
ber of wavenumber bands (or cylindrical like shells). Here, 
for k = 0.258 and A = 5 we obtain N = 25 shells. In 
real space, the total information can be recovered by summing 
over the inverse Fourier transform of each shell filtered con- 
tribution. 

The triple shell transfer occurring between the shell filtered 
quantities can be computed as, 
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and represent the basic information available to us for analy- 
sis. Knowing S > '® allows us to compute all other relevant 
nonlinear transfer quantities. By summing over all possible 
shells q we can obtain the shell-to-shell transfer (V K F \ im- 
plicitly defined below and analyzed previously 19[]) and the 
nonlinear transfer spectra by summing furthermore over p, 



P,Q 



(7) 



Numerically, when summing the transfer (jTj over k, which is 
equivalent to integrating the nonlinear transfer over the entire 
space, we obtain zero (comparable to machine precision). 

The spectral density contributions entering in the free en- 
ergy balance equation, for perpendicular characteristic scales 
k± = ItK, are presented in Fig.Q] It is interesting to note that 
while the spectral density £q is found to be zero, the £ k term, 
although it integrates to zero globally, contributes to the over- 
all linear term spectral form for time saturated states. This 
is important as the nonlinear transfer spectral density T fc is 
balanced by the sum of all the linear terms. The subsequent 
nonlinear transfers between scales can be seen as taking place 
under the constrained of a given transfer spectra. From this 
picture, the presence of the dissipation term at all scales is 
obvious. 

Locality functions. — The locality functions are defined 
from the triple transfers as a way to measure the non-locality 
degree of the triads which contribute to the energy scale flux. 
The flux through a scale (here, shell boundaries k c ) is defined 
by partial summing the transfer spectra T K , 
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In Fig. |2] we show the free energy flux across the perpendic- 
ular shell wavenumbers k±. Since the source term contribu- 
tion Q k is spread over a large interval, the flux across a scale 
k± builds up slowly to its cascade saturated value. Moreover, 
since the dissipation range is quite wide and permeates into 
the injection range, a true inertial range flux value can not be 
identified as the plateau on the flux, Fig [2] In fact, the scale 
flux plateau level is given by C + , representing the sum of the 
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FIG. 1. The free energy rhs terms spectra normalized by the total dis- 
sipation rate D. The vertical dashed lines represent the shell bound- 
aries. 
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FIG. 2. The free energy flux across the shell boundaries normalized 
by the total dissipation rate T>. The insert picture depicts the shell- 
to-shell transfer for this run; for details see Ref. [9|]. 



positive part of the linear contribution Q k + C k + T> k , here the 
first 10 shells. The C + /T> ratio clearly shows that only a frac- 
tion (57%) of the energy injected into the system contributes 
to the nonlinear cascade. 

Knowing the flux, the infrared (IR) locality function is de- 
fined by taking a probe wavenumber boundary k p , so that 
k p < k c , (sums considered in operatorial sense, as abbrevi- 
ation, due to space limitations) 
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It measures the contribution to the flux through k c from triads 
of modes with at least one wavenumber less than k p . In the 
second term, the sum over shell q starts from p + 1 to avoid 
double counting. In the limit k p — > k c , we recover the flux 
across the cutoff wavenumber k c . It is customary to normal- 
ize the locality functions to the flux trough k c , in which case 
a value of one is obtained for k p = k c and less than one for 
k p /k c < 1. Although the IR functions have a clear interpre- 
tation as the ratio of energy contributed to the flux through 
scale k c coming only from larger and larger scales, it should 
be remembered that for k p /k c <C 1 the transfers can only 
take place between the most non-local triads, i.e. triads with 
one wavevectors leg much smaller compared to the other two. 
Therefore, these functions can provide information regarding 
the locality of the nonlinear interaction. 

A similar definition is made for the ultraviolet (UV) locality 
functions, k c < k p , (sums consider as for Eq. ((9), for brevity) 
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which measures the contribution to the flux through k c from 
triads of modes with at least one wavenumber greater than k p , 
therefore providing information regarding the locality makeup 
of a scale k c in relation with smaller and smaller scales. 

Looking at the plot of 11;,- (k p \k c ) /H(k c ) as a function of 
k p /k c and n uv (k p \k c ) /H(k c ) as a function of k c /k p will re- 
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FIG. 3. The IR and UV locality functions, displayed for selected 
cutoff wavenumbers identified by the shell index c. Dashed lines 
equal or proportional to different power laws of the abscissae are 
displayed for reference. 



veal information related to the locality characteristic of the 
non-linear terms. The collapse of the locality functions de- 
pendence on k p for different values of k c represents a clear 
sign of self-similarity of the nonlinear interactions, which im- 
plies a dominance of the nonlinear terms in regard to the linear 
ones. Moreover, if the mentioned collapse exhibits a slope (in 
a log-log scale), then a state of asymptotic locality can be in- 
ferred, i.e. the nonlinear interactions saturated dynamically 
to the same level no mater how large the turbulence level be- 
comes. From our simulations, Fig. [3j none of these two be- 
haviors can be clearly observed. 

Discussion and conclusions. — Theoretically, an exponent 
for the IR and UV locality functions can be determined for 
an infinitely long inertial range. Considering the v# • V/i 
form of the nonlinearity, where v# = z x V0i and z ~ 
Vx x Vy, the triple shell transfer is found to have the form 
gK,p,Q ^ i^k^Q . where the angle brackets refers 

to volume averaging. Employing similar arguments as in 
Ref. ll2lll regarding the smoothness of scale filtered quanti- 
ties, we determine a theoretical bound for the triple trans- 



fer as, S K ' P ' Q < (\h K \)(\v%\)(\\7h p \). The scaling of the 
non-adiabatic part of the ion distribution function ((\h K 

-1/6 
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was given for two dimensional GK turbulence by Ref. |8], 
from the analysis of the self-similarity statistical symmetry 
of their respective increments. Using them, the triple scale 
transfer scales as, S K,P ® < Qq pp~ '. Assuming 

a local dyadic transfer, i.e. kg = [kx/^ — kp, kx + kp], 
in the IR limit (kp < therefore &q w kx) or in the 

UV limit (kp > 2k k therefore kq k, kp) we obtain the 
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scaling, S K,P '® ~ pp"k K J '". This in turn translates as a 
(k p /k c ) ±5 / 6 scaling for the IR and UV locality functions. 

Although an asymptotic 5/6 scaling of the IR and UV lo- 
cality functions seams plausible and would indicate a more lo- 
cal interaction compared to magnetohydrodynamic turbulence 
(2/3. Refs. fllrlTll) but more non-local compared to fluid tur- 
bulence (4/3), these values can not be clearly identified from 
our simulations. First, we need to consider that the theoret- 
ical 5 /6 exponent is found in the limit of an infinite inertial 
range, an ansatz not verified in any range for GK turbulence. 
In spite of the local energy cascade Fig.[2]-(insert), due to dis- 
sipation, the interaction of a given scale with smaller ones will 
be strongly damped increasing the scaling of the UV locality 
functions. The same scale will itself be damped compared to 
the larger scales, decreasing the IR locality function exponent. 

An effective non-local IR contribution signifies a depen- 
dence of GK turbulence on the type of instability driven it, 
while a stronger local UV depicts an insensitivity of GK's on 
the small scales and therefore the type of collision mechanism 
employed. Ultimately, these effects need to be properly ac- 
counted for in any numerical simulation. 
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